      SUBROUTINE F01BRY(N,ICN,A,LICN,LENR,IDISP,IP,IQ,LENOFF,IW,IW1,
     *                  ABORT,IFAIL)
C     MARK 7 RELEASE. NAG COPYRIGHT 1978
C     MARK 11.5(F77) REVISED. (SEPT 1985.)
C     DERIVED FROM HARWELL LIBRARY ROUTINE MC23A
C
C     PERMUTES THE ROWS AND COLUMNS OF A SPARSE MATRIX SO THAT IT IS
C     IN BLOCK LOWER TRIANGULAR FORM
C
C     .. Scalar Arguments ..
      INTEGER           IFAIL, LICN, N
      LOGICAL           ABORT
C     .. Array Arguments ..
      DOUBLE PRECISION  A(LICN)
      INTEGER           ICN(LICN), IDISP(10), IP(N), IQ(N), IW(N,5),
     *                  IW1(N,2), LENOFF(N), LENR(N)
C     .. Local Scalars ..
      INTEGER           I, I1, I2, IBEG, IBLOCK, IEND, II, ILEND, INEW,
     *                  IOLD, IROWB, IROWE, ISAVE, J, JJ, JNEW, JNPOS,
     *                  JOLD, K, LARGE, LENI, NERR, NUM, NUMNZ, NZ
C     .. Local Arrays ..
      CHARACTER*56      REC(3)
C     .. External Subroutines ..
      EXTERNAL          F01BRV, F01BRX, X04AAF, X04BAF
C     .. Intrinsic Functions ..
      INTRINSIC         MAX, MIN, MOD
C     .. Executable Statements ..
      ISAVE = IFAIL
      IFAIL = 0
C
C     SET UP POINTERS IW(.,1) TO THE BEGINNING OF THE ROWS AND SET
C     LENOFF EQUAL TO LENR.
      IW1(1,1) = 1
      LENOFF(1) = LENR(1)
C     INITIALIZE NUM AND LARGE TO PREVENT ASSIGNMENT OF UNSET
C     VARIABLES AT END OF ROUTINE IN THE EVENT OF A FAILURE.
      NUM = 0
      LARGE = 0
      IF (N.EQ.1) GO TO 40
      DO 20 I = 2, N
         LENOFF(I) = LENR(I)
         IW1(I,1) = IW1(I-1,1) + LENR(I-1)
   20 CONTINUE
C     IDISP(1) POINTS TO THE FIRST POSITION IN A/ICN AFTER THE
C     OFF-DIAGONAL BLOCKS AND UNTREATED ROWS.
   40 IDISP(1) = IW1(N,1) + LENR(N)
C
C     FIND ROW PERMUTATION IP TO MAKE DIAGONAL ZERO-FREE.
      CALL F01BRX(N,ICN,LICN,IW1,LENR,IP,NUMNZ,IW)
C
C     POSSIBLE ERROR RETURN FOR STRUCTURALLY SINGULAR MATRICES.
      IF (NUMNZ.NE.N .AND. ABORT) GO TO 340
C
C     IW1(.,2) AND LENR ARE PERMUTATIONS OF IW1(.,1) AND
C     LENR/LENOFF SUITABLE FOR ENTRY
C     TO F01BRV SINCE MATRIX WITH THESE ROW POINTER AND LENGTH
C     ARRAYS HAS MAXIMUM NUMBER OF NON-ZEROS ON THE DIAGONAL.
      DO 60 II = 1, N
         I = IP(II)
         IW1(II,2) = IW1(I,1)
         LENR(II) = LENOFF(I)
   60 CONTINUE
C
C     FIND SYMMETRIC PERMUTATION IQ TO BLOCK LOWER TRIANGULAR FORM.
      CALL F01BRV(N,ICN,LICN,IW1(1,2),LENR,IQ,IW(1,4),NUM,IW)
C
      IF (NUM.NE.1) GO TO 120
C
C     ACTION TAKEN IF MATRIX IS IRREDUCIBLE.
C     WHOLE MATRIX IS JUST MOVED TO THE END OF THE STORAGE.
      DO 80 I = 1, N
         LENR(I) = LENOFF(I)
         IP(I) = I
         IQ(I) = I
   80 CONTINUE
      LENOFF(1) = -1
C     IDISP(1) IS THE FIRST POSITION AFTER THE LAST ELEMENT IN THE
C     OFF-DIAGONAL BLOCKS AND UNTREATED ROWS.
      NZ = IDISP(1) - 1
      IDISP(1) = 1
C     IDISP(2) IS THE POSITION IN A/ICN OF THE FIRST ELEMENT IN THE
C     DIAGONAL BLOCKS.
      IDISP(2) = LICN - NZ + 1
      LARGE = N
      IF (NZ.EQ.LICN) GO TO 400
      DO 100 K = 1, NZ
         J = NZ - K + 1
         JJ = LICN - K + 1
         A(JJ) = A(J)
         ICN(JJ) = ICN(J)
  100 CONTINUE
      GO TO 400
C
C     DATA STRUCTURE REORDERED.
C
C     FORM COMPOSITE ROW PERMUTATION ... IP(I) = IP(IQ(I)).
  120 DO 140 II = 1, N
         I = IQ(II)
         IW(II,1) = IP(I)
  140 CONTINUE
      DO 160 I = 1, N
         IP(I) = IW(I,1)
  160 CONTINUE
C
C     RUN THROUGH BLOCKS IN REVERSE ORDER SEPARATING DIAGONAL
C     BLOCKS WHICH ARE MOVED TO THE END OF THE STORAGE.  ELEMENTS IN
C     OFF-DIAGONAL BLOCKS ARE LEFT IN PLACE UNLESS A COMPRESS IS
C     NECESSARY.
C
C     IBEG INDICATES THE LOWEST VALUE OF J FOR WHICH ICN(J) HAS
C     BEEN SET TO ZERO WHEN ELEMENT IN POSITION J WAS MOVED TO THE
C     DIAGONAL BLOCK PART OF STORAGE.
      IBEG = LICN + 1
C     IEND IS THE POSITION OF THE FIRST ELEMENT OF THOSE TREATED
C     ROWS WHICH ARE IN DIAGONAL BLOCKS.
      IEND = LICN + 1
C     LARGE IS THE DIMENSION OF THE LARGEST BLOCK ENCOUNTERED SO
C     FAR.
      LARGE = 0
C
C     NUM IS THE NUMBER OF DIAGONAL BLOCKS.
      DO 300 K = 1, NUM
         IBLOCK = NUM - K + 1
C        I1 IS FIRST ROW (IN PERMUTED FORM) OF BLOCK IBLOCK.
C        I2 IS LAST ROW (IN PERMUTED FORM) OF BLOCK IBLOCK.
         I1 = IW(IBLOCK,4)
         I2 = N
         IF (K.NE.1) I2 = IW(IBLOCK+1,4) - 1
         LARGE = MAX(LARGE,I2-I1+1)
C        GO THROUGH THE ROWS OF BLOCK IBLOCK IN THE REVERSE ORDER.
         DO 280 II = I1, I2
            INEW = I2 - II + I1
C           WE NOW DEAL WITH ROW INEW IN PERMUTED FORM (ROW IOLD IN
C           ORIGINAL MATRIX).
            IOLD = IP(INEW)
C           IF THERE IS SPACE TO MOVE UP DIAGONAL BLOCK PORTION
C           OF ROW GO TO ...
            IF (IEND-IDISP(1).GE.LENOFF(IOLD)) GO TO 220
C
C           IN-LINE COMPRESS.
C           MOVES SEPARATED OFF-DIAGONAL ELEMENTS AND UNTREATED ROWS TO
C           FRONT OF STORAGE.
            JNPOS = IBEG
            ILEND = IDISP(1) - 1
            IF (ILEND.LT.IBEG) GO TO 360
            DO 180 J = IBEG, ILEND
               IF (ICN(J).EQ.0) GO TO 180
               ICN(JNPOS) = ICN(J)
               A(JNPOS) = A(J)
               JNPOS = JNPOS + 1
  180       CONTINUE
            IDISP(1) = JNPOS
            IF (IEND-JNPOS.LT.LENOFF(IOLD)) GO TO 360
            IBEG = LICN + 1
C           RESET POINTERS TO THE BEGINNING OF THE ROWS.
            DO 200 I = 2, N
               IW1(I,1) = IW1(I-1,1) + LENOFF(I-1)
  200       CONTINUE
C
C           ROW IOLD IS NOW SPLIT INTO DIAG. AND OFF-DIAG. PARTS.
  220       IROWB = IW1(IOLD,1)
            LENI = 0
            IROWE = IROWB + LENOFF(IOLD) - 1
C           BACKWARD SCAN OF WHOLE OF ROW IOLD (IN ORIGINAL MATRIX).
            IF (IROWE.LT.IROWB) GO TO 260
            DO 240 JJ = IROWB, IROWE
               J = IROWE - JJ + IROWB
               JOLD = ICN(J)
C              IW(.,2) HOLDS THE INVERSE PERMUTATION TO IQ.
C              ..... IT WAS SET TO THIS IN F01BRV.
               JNEW = IW(JOLD,2)
C              IF (JNEW.LT.I1) THEN ....
C              ELEMENT IS IN OFF-DIAGONAL BLOCK AND SO IS LEFT IN SITU.
               IF (JNEW.LT.I1) GO TO 240
C              ELEMENT IS IN DIAGONAL BLOCK AND IS MOVED TO THE
C              END OF THE STORAGE.
               IEND = IEND - 1
               A(IEND) = A(J)
               ICN(IEND) = JNEW
               IBEG = MIN(IBEG,J)
               ICN(J) = 0
               LENI = LENI + 1
  240       CONTINUE
            LENOFF(IOLD) = LENOFF(IOLD) - LENI
  260       LENR(INEW) = LENI
  280    CONTINUE
         IP(I2) = -IP(I2)
  300 CONTINUE
C     RESETS IP(N) TO POSITIVE VALUE.
      IP(N) = -IP(N)
C     IDISP(2) IS POSITION OF FIRST ELEMENT IN DIAGONAL BLOCKS.
      IDISP(2) = IEND
C
C     THIS COMPRESS IS USED TO MOVE ALL OFF-DIAGONAL ELEMENTS TO
C     THE FRONT OF THE STORAGE.
      IF (IBEG.GT.LICN) GO TO 400
      JNPOS = IBEG
      ILEND = IDISP(1) - 1
      DO 320 J = IBEG, ILEND
         IF (ICN(J).EQ.0) GO TO 320
         ICN(JNPOS) = ICN(J)
         A(JNPOS) = A(J)
         JNPOS = JNPOS + 1
  320 CONTINUE
C     IDISP(1) IS FIRST POSITION AFTER LAST ELEMENT OF OFF-DIAGONAL
C     BLOCKS.
      IDISP(1) = JNPOS
      GO TO 400
C
C     ERROR RETURN
C     IFAIL = 1 - MATRIX IS STRUCTURALLY SINGULAR, RANK = NUMNZ
C     IFAIL = 7 - LICN NOT BIG ENOUGH, INCREASE BY N
C
  340 IFAIL = 1
      GO TO 380
  360 IFAIL = 7
  380 CONTINUE
C     ** CODE FOR OUTPUT OF ERROR MESSAGES *************************
      IF (MOD(ISAVE/10,10).EQ.0) GO TO 400
      CALL X04AAF(0,NERR)
      IF (IFAIL.EQ.1) THEN
         WRITE (REC,FMT=99999) NUMNZ
         CALL X04BAF(NERR,REC(1))
         CALL X04BAF(NERR,REC(2))
         CALL X04BAF(NERR,REC(3))
      END IF
      IF (IFAIL.EQ.7) THEN
         WRITE (REC,FMT=99998) N
         CALL X04BAF(NERR,REC(1))
         CALL X04BAF(NERR,REC(2))
      END IF
C     ** END OF CODE FOR OUTPUT OF ERROR MESSAGES ******************
  400 IDISP(8) = NUMNZ
      IDISP(9) = NUM
      IDISP(10) = LARGE
      RETURN
C
99999 FORMAT (/' MATRIX IS STRUCTURALLY SINGULAR - DECOMPOSITION ABORT',
     *  'ED',/'  RANK = ',I6)
99998 FORMAT (/' LICN NOT BIG ENOUGH FOR PERMUTATION - INCREASE BY',I6)
      END
